*****************************************************************
* Replication directory for                                   ***
* Prime locations                                             ***
* by Gabriel M. Ahlfeldt, Thilo N.H. Albers, Kristian Behrens ***
* Published in American Economic Review: Insights             ***
*****************************************************************
* 01/2025
* Stata
version 17.0

* This do files generate BDD plots of the change in employment density at the prime location boundary by p-value
* The analysis will guide the choice of the preferred p-value

* Load data
	use "$temp/USMETROS_PL_distances.dta", clear 
			egen metro_emp = sum(emp), by(cbsafp)

* Prepare variables names
	gen metro_id=cbsafp
	
* RESTRICT THE DATA TO DEVELOPABLE AREA
 	drop if devel  == 0

* calculate employment density
	gen cellarea=area_g
	gen empl_density=emp/cellarea/10^4
	gen lempl_density = ln(empl_density+1)
	gen empl_density_TS = cell_TS_emp/cellarea/10^4
	gen lempl_density_TS  = ln(empl_density_TS +1)	
		
* Main loop over p-values begins ***********************************************
	foreach pvalue in  $plist { 
		BDD "lempl_density_TS empl_density_TS lempl_density empl_density" "dist_`pvalue'"	"`pvalue'" 0.100 21.00 1
	}
	
* End of main loop over p-values ***********************************************		  

* Plot results
cd "$temp"
	graph combine	"_FIG_lempl_density_p98_5.gph"  ///
					"_FIG_lempl_density_p99_0.gph"  ///
					"_FIG_lempl_density_p99_5.gph"  ///
					"_FIG_lempl_density_p99_9.gph"  ///
					, graphregion(color(white) lwidth(thick)) scale(0.75) ycommon ///
					l1("Log employment density (in 10,000 Jobs/km{superscript:2})") b1("Distance from PL border (in km)") xsize(10) ysize(5) 
* Export Appendix Figure B.1.1, panel a	
	capture mkdir "$figures_App/US-CBSAs"
	graph export "$figures_App/US-CBSAs/FIG_B1_1a_PctCutoff-log.pdf",  replace	
cd "$root"
		
* Restrict to Top-10 cities 
	preserve
	keep metro_id metro_emp
	duplicates drop
	gsort -metro_emp
	local T = metro_emp[10]
	restore
	keep if metro_emp >= `T'
		
* Main loop over p-values begins ***********************************************
	foreach pvalue in  $plist { 
		BDD "lempl_density_TS empl_density_TS lempl_density empl_density" "dist_`pvalue'"	"`pvalue'" 0.100 21.00 1
	}
	
* End of main loop over p-values ***********************************************		

* Plot results

cd "$temp"
	graph combine	"_FIG_lempl_density_p98_5"  ///
					"_FIG_lempl_density_p99_0"  ///
					"_FIG_lempl_density_p99_5"  ///
					"_FIG_lempl_density_p99_9"  ///
					, graphregion(color(white) lwidth(thick)) scale(0.75) ycommon ///
					l1("Log employment density (in 10,000 Jobs/km{superscript:2})") b1("Distance from PL border (in km)") xsize(10) ysize(5)  
* Export Appendix Figure B.1.1, panel b	
	graph export "$figures_App/US-CBSAs/FIG_B1_1b_PctCutoff-log-top10.pdf",  replace	

cd "$root"
* Script ends
